pacman::p_load(tidyverse, ggrepel, FactoMineR, magrittr)
setwd("D:/project/R-use/yunan/script/20201116")


da <- function(da) {
  (da / (da %>%
    colSums() %>%
    as.numeric() %>%
    rep(each = nrow(da)) %>%
    matrix(ncol = ncol(da), nrow = nrow(da)))) %>%
    t() %>%
    return()
}

da <- "data.xls" %>%
  read.delim(header = T, check.names = F, row.names = 1, sep = "\t") %>%
  da()
#载入分组数据
groups <- "group2.txt" %>% read.delim(header = T, sep = "\t")

pca <- da %>% PCA(scale.unit = F, graph = F)
pca %$% write.csv(ind$coord, file = "pca.xls")

group <- NULL

for (i in seq_along(groups$Sample)) {
  group[i] <- groups$Group[rownames(pca$ind$coord)[i] %>% grep(groups$Sample)]
}


if (ncol(groups) == 2) {
  pcadata <- data.frame(rownames(pca$ind$coord), pca$ind$coord[, 1], pca$ind$coord[, 2], group) %>%
    set_colnames(c("sample", "PC1", "PC2", "group"))
}

ggplot(pcadata, aes(PC1, PC2)) +
  geom_point(aes(fill = group, color = group), size = 4) +
  labs(title = "PCA Plot",
       x = (floor(pca$eig[1, 2] * 100) / 100) %>% paste0("PC1 ( ", ., "%", " )"),
       y = (floor(pca$eig[2, 2] * 100) / 100) %>% paste0("PC2 ( ", ., "%", " )")) +
  theme(text = element_text(family = "Times", size = 10)) +
  theme(panel.background = element_rect(fill = 'white', colour = 'black'),
        axis.title.x = element_text(colour = "black", size = 10),
        axis.title.y = element_text(colour = "black", size = 10),
        axis.text = element_text(colour = "black", size = 10),
        legend.title = element_text(colour = "black", size = 10, face = "bold"),
        legend.key = element_blank(),
        plot.title = element_text(size = 10, colour = "black", hjust = 0.5, face = "bold")) +
  theme(legend.position = c(0.9, 0.9)) +
  theme(legend.title = element_blank())+
  geom_vline(aes(xintercept = 0), linetype = "dotted") +
  geom_hline(aes(yintercept = 0), linetype = "dotted")


